## This file provides the code to produce the tables and figures in:
## Sharan Grewal, "Military Defection during Localized Protests: The Case of
## Tataouine," International Studies Quarterly, Forthcoming.


## Load data and packages

library(stargazer)
library(ggplot2)
library(scales)
library(sampleSelection)
library(MASS)
data <- read.csv("data.csv")


## Table 1: Descriptive Statistics
table(data$rank)/72
table(data$army)/72
table(data$navy)/72
table(data$airforce)/72
table(data$milsec)/72
table(data$shared)/72
table(data$region)/72
table(data$train_US)/72
table(data$train_France)/72
table(data$train_other)/72



## Table 2:
one <- lm(defect~interior+NSC+noparty, data=data)
two <- lm(defect~interior+NSC+noparty+rank+army+retrev+train_US+essebsi, data=data)
three <- lm(defect2~interior+NSC+noparty+rank+army+retrev+train_US+essebsi, data=data)
stargazer(one, two, three,
          dep.var.labels=c("Exc. DK (1-4)", "Inc. DK (1-5)"),
          covariate.labels=c("Interior", "NSC", "Apolitical", "Rank", "Army", "Ret. Year", "U.S. Training", "Essebsi"))


## Figure 1: Composition
data.summary <- predict(one, newdata=data.frame(interior=c(0,1), NSC=rep(mean(data$NSC, na.rm=T),2), noparty=rep(mean(data$noparty, na.rm=T), 2)), interval="confidence", level=0.9)
data.summary <- data.frame(treatment=c("Coast (N=43)", "Interior (N=29)"), data.summary)
ggplot(data.summary, 
       aes(x=data.summary$treatment, y=data.summary$fit)) +  
  geom_bar(position=position_dodge(width=0.9), stat="identity", fill=c("darkgrey", "darkblue")) + 
  geom_errorbar(aes(ymin=data.summary$lwr, 
                    ymax=data.summary$upr, width=0.1),
                position=position_dodge(width=0.9)) +
  ggtitle("Effect of Composition on Defection") +
  theme_bw() + 
  theme(panel.grid.major = element_blank()) + 
  theme(plot.title = element_text(hjust = 0.5)) +
  scale_y_continuous(limits=c(1, 5),breaks=seq(1, 4, by=1),oob=rescale_none) +
  xlab("") +
  ylab("Support for Defection (1-4)") +
  geom_text(data=data.summary, 
            aes(x=data.summary$treatment, y=data.summary$fit,
                label=round(data.summary$fit,2)),
            vjust=-5.5,size=5,position=position_dodge(0.9)) +
  theme(text = element_text(size=16)) +
  geom_text(aes(x=1.5, y=4.8, label="p=0.003"), size=5)


## Figure 2: Corporate interests
data.summary <- predict(one, newdata=data.frame(interior=rep(mean(data$interior, na.rm=T),3), NSC=c(3,4,5), noparty=rep(mean(data$noparty, na.rm=T), 3)), interval="confidence", level=0.9)
data.summary <- data.frame(treatment=c("Neutral (N=6)", "Agree (N=31)", "Strongly Agree (N=33)"), data.summary)
data.summary$treatment <- factor(data.summary$treatment, levels=c("Neutral (N=6)", "Agree (N=31)", "Strongly Agree (N=33)"))
data.summary$group <- 1
ggplot(data.summary, 
       aes(x=data.summary$treatment, y=data.summary$fit)) +  
  geom_point() + 
  geom_line(group=1) +
  geom_errorbar(aes(ymin=data.summary$lwr, 
                    ymax=data.summary$upr, width=0.1),
                position=position_dodge(width=0.9)) +
  ggtitle("Effect of Corporate Interests on Defection") +
  theme_bw() + 
  theme(panel.grid.major = element_blank()) + 
  theme(plot.title = element_text(hjust = 0.5)) +
  scale_y_continuous(limits=c(1, 4.5),breaks=seq(1, 4, by=1),oob=rescale_none) +
  xlab("\nMore Military Reps in NSC") +
  ylab("Support for Defection (1-4)") +
  geom_text(data=data.summary, 
            aes(x=data.summary$treatment, y=data.summary$fit,
                label=round(data.summary$fit,2)),
            vjust=c(-5.5,-4,-4.5),size=5,position=position_dodge(0.9)) +
  theme(text = element_text(size=16)) +
  geom_text(aes(x=2, y=4.3, label="p=0.03"), size=5)



## Table 3: Heckman Selection Models

summary(selection(notDK~interior+NSC+noparty+train_US+rank+army+ryear+essebsi, defect~interior+NSC+noparty+train_US+rank+ryear+essebsi, data=data))
one <- selection(notDK~interior+NSC+noparty+train_US+army+rank+ryear+essebsi, defect~interior+NSC+noparty+train_US+rank+ryear+essebsi, data=data)
two <- selection(notDK~interior+NSC+noparty+train_US+army+rank+ryear+essebsi, defect~interior+NSC+noparty+army+rank+ryear+essebsi, data=data)
three <- selection(notDK~interior+NSC+noparty+train_US+army+rank+ryear+essebsi, defect~interior+NSC+train_US+army+rank+ryear+essebsi, data=data)
stargazer(one,two,three,
          dep.var.labels=c("Exc. DK (1-4)", "Exc. DK (1-4)", "Inc. DK (1-5)"),
          covariate.labels=c("Interior", "NSC", "Apolitical", "U.S. Training", "Army","Rank", "Ret. Year", "Essebsi"))




## Appendix D: Ordered Logit

one <- polr(lm(factor(defect)~interior+NSC+noparty, data=data), method="logistic")
two <- polr(lm(factor(defect)~interior+NSC+noparty+rank+army+train_US+essebsi, data=data), method="logistic")
three <- polr(lm(factor(defect2)~interior+NSC+noparty+rank+army+retrev+train_US+essebsi, data=data), method="logistic")
stargazer(one, two, three,
          dep.var.labels=c("Exc. DK (1-4)", "Exc. DK (1-4)", "Inc. DK (1-5)"),
          covariate.labels=c("Interior", "NSC", "Apolitical", "Rank", "Army", "Ret. Year", "U.S. Training", "Essebsi"))


## Questions? Contact Sharan Grewal at ssgrewal@princeton.edu


